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The constraints of a realistic magnetic interaction limit the possibility of realizing degenerate 
vertex models with artificial spin ice in two dimensions. Here we present and study a novel vertex- 
frustrated geometry which exhibits a residual entropy density. While interesting in its own regard as 
a new form of frustration in a vertex system, our lattice opens new pathways for implementing de- 
generacy in artificial magnetic nanoarrays, a task so far elusive. We demonstrate that its degenerate 
ground state can be mapped to the disordered phase of the exactly solvable _F-model. Monte Carlo 
simulations extend our study to thermal states, and, significantly, find algebraic correlations for 
frustrated vertices, magnetic charges, and spins in the ground state, pointing to a critical phase. 



Artificial Spin Ice (ASI) has raised considerable inter- 
est for its technological potentials, and as a tailorable 
medium to investigate novel collective phenomena in a 
materials- by-design approach [l~18|. ASI was inspired 
by the so-called spin ice compounds [l^, [1^, a special 
class of pyrochlore ferromagnets in which the magnetic 
interactions resemble the frustrated proton ordering in 
water ice [ill - EsI l . Both water ice and spin ice retain a fi- 
nite entropy density even at very low temperatures. The 
local spin arrangement in spin ice is dictated by rules 
similar to the Bernal-Fowler ice rules [2lj. An interest- 
ing consequence of this local ordering is the appearance 
of algebraically decaying spin correlations (23 | , similar to 
those in two-dimensional (2D) lattice systems with hard- 
core constraints such as dimcr or vertex models [25 1. 

The original artificial spin ice presented by Wang et al. 
consists of magnetically interacting elongated permalloy 
nanoislands arranged as links of a square lattice P, SIJ. 
Within the nearest-neighbor approximation, the square 
ice can be mapped to a non-degenerate six-vertex model: 
the anisotropic nature of the magnetic interaction in a 
2D arrangement lifts the degeneracy of the six distinct 
2-in-2-out vertices and the ground state is described by 
the ordered phase of the F-model |25| . Ways to overcome 
this limit and induce degeneracy have been proposed (26j 
yet they present challenges in nanofabrication. 

So far, kagome ice, in which the islands are arranged 
along the edges of a honeycomb lattice 12 - 13] , is the only 
ASI of extensive degeneracy, inheriting residual entropy 
from its six low-energy 2-in-l-out or l-in-2-out degener- 
ate vertices (pseudo-ice rule). However this pseudo-ice 
regime is non-critical, with exponentially-decaying spin 
correlation [3l|. Recently Morrison and coworkers have 
proposed "vertex-frustrated lattices" for ASI. The resid- 
ual entropy in these systems emerges from the inability 
of allocating all the vertices in their lowest energy con- 
figuration, rather than from a degeneracy at the vertex 
level 
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Here we present the first comprehensive study of a 
particularly interesting vertex-frustrated lattice depicted 
in Fig. [TJa). This so-called shakti lattice is maximally 
vertex- frustrated [isj : on each minimal loop, at least one 
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FIG. 1: (Color online) (a) The lattice geometry of the shakti 
spin ice. The spin-ice state and the corresponding defect- 
vertex configuration in a typical disordered ground state are 
shown in panels (b) and (c), respectively. The defect type-IF 
vertices are indicated by circles. The defect configuration is 
further mapped to an emergent 6- vertex model in (d). 



vertex is not in its lowest energy configuration. Remark- 
ably, we find that its residual entropy can be computed 
exactly by mapping the T = manifold of the shakti 
spin-ice into a disordered state of the non-degenerate F- 
model. The degeneracy implied by topological frustra- 
tion can be subsumed as a fictitious temperature for a 
non-degenerate vertex model. We also show numerically 
that its correlations are algebraic, indicating a critical 
ground state. 

The shakti lattice can be derived from the square lat- 
tice [l| by alternatively placing an additional vertical 
or horizontal permalloy island in each square plaquette. 
The resulting lattice consists of vertices with both coordi- 
nation numbers z = 4 and 3, as in Fig.[2j The perpendic- 
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FIG. 2; (Color online) The seven different types of vertices 
in artificial spin ice shown in Figs. [Ha) and (b) . Also shown 
is the magnetostatic energy of the vertices expressed in terms 
of three parameters a, a', and /?. The multiplicities of the 
various vertex types are: qi — 2, qu = 4, qm — 8, and giv ~ 2 
for 2 = 4 vertices, and qy = 2, qii' = 4, and qmi = 2 for z = 3 
vertices. 



ular geometry of the islands around the z = 3 sites breaks 
the six-fold degeneracy of the 2-in-l-out/l-in-2-out ver- 
tices. A new type of frustration arises from the fact that 
not all z = 3 vertices can assume the minimum-energy 
configurations simultaneously. 

In this paper we limit ourselves to a vertex-model of 
the shakti lattice, which has proved useful in describing 
experimental results in ASI [3, S Q • Energetically, there 
are four distinct vertex types, labeled by numerals I-IV, 
for the z = 4 sites, while the three types of z = 3 ver- 
tices are labeled as type F, IF, and IIF. The dominant 
nearest-neighbor interaction energy can be parameter- 
ized by a and /?, the magnetostatic energy between two 
islands in perpendicular and parallel geometries, respec- 
tively, at the z — A vertex, and a' represents the energy 
between a pair of perpendicular island at the z = 3 ver- 
tex. In terms of these parameters, the various vertex 
energies are summarized in Fig. [5] 

Since the energies a > /3 in 2D, contrary to the py- 
rochlore spin ice, the degeneracy between the two types 
(I and II) of 2-in-2-out vertices is lifted. In the case of 
square ASI, a long-range spin order consisting of a stag- 
gered arrangement of the energetically favorable type-I 
vertices emerges as the ground state. 

Approximating each permalloy island as a point-like 
dipole, simple calculation gives a = 24V2, (3 = 16, and 
a' = X (384/25-\/5) measured in units of = fiQfJ-^ /4:na^ , 
where ^ and xfj. are the effective dipole moments of the 
blue and red spins, respectively, and a is the lattice con- 
stant. The geometry and size of the red spins thus deter- 
mine which vertex type (F or IF) has lower energy. As 
can be easily checked, an energetically favorable type-IF 
vertex gives rise to a long-range spin order in the shakti 
lattice. In the following, we mainly concern ourselves 
with the more interesting case in which ep < eip . 

Ideally the ground state of the spin ice is reached when 
the individual vertex assumes either type-I or F con- 
figurations depending on its coordination number [l5| . 
However, the special geometry of the lattice precludes 
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FIG. 3: (Color online) Monte Carlo simulation of the shakti 
spin ice with A'^s = 1000 spins. Parameters a = a — 1.41/3 
and periodic boundary conditions are used in the simulations, 
(a) and (b) show the fraction of various vertex types as a func- 
tion of temperature T measured in units of — fj,ofj,^ /Aiva'^ . 
(c) shows the temperature dependence of entropy density s 
and specific-heat c. The entropy curve s{T) is obtained by 
integrating the specific-heat c{T)/T. The averaged magni- 
tude of magnetic charge (in natural units) at the z = 3 and 4 
vertices is shown in (d). 



such an ideal arrangement [15|. Assuming a' ^ a, the 
type-IF is the next favorable vertex energetically. Con- 
sequently, the ground state is expected to contain only 
the two minimum- energy vertices (I and F) and the next 
lowest-energy type-IF vertices. 

To investigate the thermodynamic behavior of the 
shakti spin ice as a vertex model, we first perform Monte 
Carlo simulations. The results shown in Figs. [3] are ob- 
tained assuming nearest-neighbor interaction with pa- 
rameters a = a' = 1.41/3. At high temperatures the 
system is in an uncorrelated paramagnetic phase as the 
population of various vertex types reaches an equilibrium 
according to their respective multiplicities. On the other 
hand, as temperature tends to zero, the magnet is indeed 
dominated by vertices of types I, F, and IF. The crossover 
between these two regimes is marked by a broad peak in 
the specific-heat [Fig. Eljc)]. The presence of ice rules 
at low T is confirmed by the temperature dependence of 
averaged magnetic charges at individual vertex shown in 
Fig.[3ljd). The 2-in-2-out rule indicates a vanishing mag- 
netic charge Q = at the z = A vertices while the 1-in- 
2-out/2-in-l-out rule implies a magnetic charge Q = ±1 
at the z ~ 3 vertices. 

By fixing the high-T entropy density to fcBln2, we 
find a residual entropy sq « 0.1178 fcs per spin from 
the simulations. This nonzero entropy density implies 
an extensive degeneracy of spin-ice states at low temper- 
atures. Below the crossover, we notice that the only ver- 



tex types with nonzero fractional populations are nj = 1 , 
ni' = nil' = 1/2: while all z = 4 sites are minimum- 
energy type-I vertices, only half of the z = 3 vertices are 
in the ground state. The other half are frustrated type- 
IP vertices. We can thus use locations of the frustrated 
vertices to characterize the spin-ice states; these unhappy 
vertices will be referred to as defects for simplicity. For 
example, Figs. [U^b) and (c) show a generic disordered 
spin-ice state and the corresponding defect configuration, 
respectively. It should be noted that this mapping is at 
least 2-to-l: each spin-ice state and its time- reversal part- 
ner arc mapped to the same defect configuration. For 
square plaquette with defects sitting at the two ends of 
the center spin, i.e. type-6 plaquette in Fig. HI there is 
an additional Z2 degree of freedom associated with the 
direction of the center (red) spin. 

To understand the structure of the degenerate spin- 
ice manifold, we first note that each plaquette has at 
least two defects in any spin-ice state; this can be easily 
demonstrated by enumerating all possible magnetic con- 
figurations in a plaquette [15[. As a result, any configu- 
ration in which each plaquette has exactly two defects is 
a ground state, consistent with the numerical result that 
nw = l/2asT^-0. This two-defects constraint is sim- 
ilar to the 2-in-2-out ice rule in square ice. Fig. |4] shows 
the six different two-defects configurations in a plaque- 
tte. Since there are two types of plaquette in our lattice 
with different orientations of the center island, those with 
vertical (horizontal) center spin are referred to in the fol- 
lowing as type-A (B) plaquettes. 

Noting that the centers of plaquettes in our model form 
a square lattice, it is not surprising that each defect con- 
figuration can be mapped to a 6-vcrtex state. To es- 
tablish this mapping, we draw an arrow from type-A to 
its neighboring typc-B plaquette if their common edge 
contains a defect vertex. On the other hand, an edge 
without defect vertex corresponds to an arrow pointing 
from type-B to type-A plaquette. This establishes an 
one-to-one mapping between the two-defects plaquettes 
and the 2-in-2-out vertices; see Fig. SI An example of the 
mapping is also illustrated in Fig.[IJd). 

We now consider the partition function of the emer- 
gent 6-vertcx model. As discussed above. In addition to 
the global Z2 time-reversal symmetry, there are two dif- 
ferent spin states corresponding to typc-6 vertex. This 
additional doublet degeneracy comes from the fact that 
the direction of center spin can be parallel or anti-parallel 
to the orientation of the center island [Fig. UJ^f)]. Let 
denote the number of type-6 vertex in a given 6-vertex 
configuration C, the ground-state partition function of 
the shakti spin ice is simply 



(1) 



The second identity holds when periodic or cylindrical 
boundary conditions are imposed. This is because type- 
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FIG. 4: (Color online) The six possible two-defects config- 
urations in a plaquette (for both type-A and B) and their 
mapping to the 2-in-2-out vertices. 



5 and 6 vertices are sources and sinks, respectively, of 
the horizontal arrows or 'fluxes'. As the total flux is 
conserved, there must be as many sinks as sources, i.e. 
ris = rig [25| . The degeneracy of Eq. ([T]) can be estimated 
as: W = Wico X 2^n<"'=), where Wice = (|)^^°/^° is the 
degeneracy of the ideal square ice [IHj, (ng) ^ 0.189 is 
the average population of type-6 vertex in a ice state, 
and N\j is the number of plaquettes (or vertices in a 6- 
vertex state). Since Nu = Ns/5 in shakti ice, where 
Ns is the number of spins, this gives a residual entropy 
density sq = kB^nW/Ng ^ 0.1125A:b, which is close to 
the numerical result. 

The entropy density of the shakti spin-ice can be ex- 
actly computed by recasting the partition function ([T]) 
into that of the standard F- model 271: 



Z = V2^° Zf = V2' 



C i=l 



(2) 



where uji denotes the statistical weight of the type-i ver- 



tex: LUi 



UJ2 = 



W4 



-K 



and W5 = wg — 1. The 



case of shakti ice corresponds to K = ^ In 2. 

For a general F-model, the parameter K > acts as 
an inverse temperature; the square ice corresponds to the 
special case with if = or infinite T. It is known that 
a Kosterlitz-Thouless transition occurs at Kc = In 2 [2^ . 
The correlation length is finite in the ordered state for 
K > Kc, while the system remains critical with an in- 
finite correlation length for K < Kc- Our case with 
K = Kc/2 thus corresponds to a critical phase with dis- 
ordered vertices and spins. 

The residual entropy of Eq. © is S'/fc^ ~ InZ = 
(7Vn/2) In 2 -I- In Zp. The entropy density of the i^-model 
with K = Kc/2 can be exactly computed using the Bethe 
ansatz [2?! . We thus obtain a residual entropy density for 
the shakti ice: 



S 



1 



lncot{k/2) dk 



2G 



0.583122, (3) 



where G is the Catalan's constant [S^l . Since the number 
of spins Ns = 5iVn in shakti ice, Eq. ^ corresponds to 
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FIG. 5: (Color online) The correlation function for (a) defect- 
vertex, (b) magnetic charge, and (c) magnetic moment ob- 
tained from Monte Carlo simulations using the loop algo- 
rithm These correlations are measured between defects, 
charges, or spins at the vertical edges of two plaquettes sepa- 
rated by r unit cells along the diagonal direction of the lattice. 
The insets show the same data (with absolute values) in log- 
log plot. Panel (d) shows the probability distribution of loop 
lengths s. 



an entropy density per spin S/N^ ks = 0.116624, which 
is fairly close to the numerically estimated value sq = 
0.1178/cB [Fig. [Sic)]. 

We now turn to the correlation functions in the ground 
state. To efRciently navigate between different six-vertex 
states, we use the loop algorithm in which detailed bal- 
ance is always satisfied locally when constructing the 



loops [28j, |29|. We first consider the defect-defect cor- 
relation function C„(r) = ((57i(r) 5n(0)), where 6n{v) = 
n{v) — (n), and n(r) = 1 or dependeing on the pres- 
ence or absence of defect vertex at edge r. The average 
defect density (n) = 1/2 is consistent with the constraint 
that each plaquette has two defects in the ground state. 
The correlation function C„(r) along the diagonal direc- 
tion of the lattice, shown in Fig. EJa), falls off quickly 
beyond a few lattice constants. As expected for corre- 
lations in a critical state, the inset log-log plot shows 
a power- law decaying C„(r) ^ r~''" with an estimated 
exponent = 1.96986 ± 0.0046. 

As in the kagome ice case, there is a residual mag- 
netic charge Q = ±1 at the z = 3 vertices [also see 
Fig. IHd)]. The local arrangement of charges also fol- 
lows certain rules: charges of the same sign always occur 
at the two ends of the center vertical or horizontal island 
in type-1 to 4 plaquettes, while opposite charges always 
occupy the two ends of the center island in a typc-5 or 6 



plaquette. Fig. [5jb) shows a stronger charge-charge cor- 
relation Gq{y) = {Q{r) (5(0)) than that of defect vertices. 
Again, we find an algebraic correlation CQ{r) ^ r~'"5 at 
large separation with r]Q = 0.497144 ± 0.001952. 

Since the mapping to defects configuration in a plaque- 
tte and the corresponding 6- vertex state is at least 2-to- 
1, the magnetic state of a plaquette is determined by its 
vertex type and an Ising variable. For a given 6-vertex 
state C, these Ising variables arc not uniform and have 
to be determined consistently by first fixing spins at the 
boundaries. The numerically obtained spin-spin correla- 
tion function is shown in Fig. (Sfc) . Its power- law decay 
Cs{r) - r-''= has an exponent r]s = 0.502075 ± 0.002028. 

The lack of a length scale in a critical phase implies 
that loops of all lengths will be found in any given state. 
To see this we numerically evaluate the probability distri- 
bution function P{s) of loop length s; the result is shown 
in Fig. [SJ^d). Similar to loops in square or pyrochlore 
ices [30[, the loop statistics exhibits two distinct regimes: 



a power-law distribution P{s) 



-1.125 



for short loops 



and a fiat distribution due to winding loops in a finite 
system. 

In summary, we have proposed a novel artificial spin 
ice characterized by frustrated ordering of defect vertices 
in the ground state. The resulting degenerate manifold 
can be mapped to a special 6-vertex F-model and solved 
exactly. We have also solved the model numerically and 
found that in the ground state, correlation functions of 
frustrated vertices, magnetic charges, and spins exhibit 
a power-law decay consistent with a critical ground state 
predicted by the exact solvable F-model. Beside its in- 
trinsic interest as a novel form of frustration which can be 
mapped into temperature, this system can be easily fab- 
ricated as an extensively degenerate ASI, and annealed in 
its degenerate ground state 16- 1^], thus opening novel 
directions in the study of frustration induced degeneracy. 
In future work we will demonstrate the emergence of a 
smectic phase corresponding to magnetic monopolc crys- 
tallization when long range interactions are included. 
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